library(DESeq2)

A=read.csv("C:/本科生教学/大三/R/实践7/mRNA_exprSet.csv")
rownames(A)=A[,1]
A=A[,-1]
head(A)
condition <- factor(ifelse(grepl(".{13}11",colnames(A)),"normal","tumor"))
condition
colData <- data.frame(row.names=colnames(A), condition)
colData
dds <-  DESeqDataSetFromMatrix(countData = A,colData = colData,design = ~ condition)
dds1=DESeq(dds)   # DEseq2差异分析
res=results(dds1)   #提取结果
B=as.data.frame(res)   #将差异分析结果转换为数据框
write.csv(B,file = "实践7DESeq.csv")   #生成结果文件
B=subset(B, !is.na(padj) & !is.na(log2FoldChange))   #去除有缺失值不便后续分析的基因


library(ggplot2)
library(ggrepel)
B$group=as.factor(ifelse(B$padj<0.05 & B$log2FoldChange>log2(2),"up",
                               ifelse(B$padj<0.05 & B$log2FoldChange<log2(1/2),"down","NS")))
volcanoplot=ggplot(B, aes(x=log2FoldChange, y=-log10(padj)))+
  geom_point(aes(color = B$group), alpha=0.7, size=2)+    #点透明度和大小
  scale_color_manual(values=c("blue","grey","red"))+
  labs(title="volcanoplot",
       x="log2FoldChange", 
       y="-log10(padj)")+
  geom_hline(yintercept = -log10(0.05),linetype = 'dashed', size = 0.2)+        # 在图上添加虚线
  geom_vline(xintercept=c(-log2(2),log2(2)),linetype = 'dashed', size = 0.2)# 在图上添加虚线
ggsave(volcanoplot,filename = "volcanoplot实践7.jpg")
volcanoplot
dev.off()
